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Abstract 

The initial boundary value problem describing the evolution of unsteady 
linearized perturbations of a steady, uniform subsonic flow is analyzed. The 
eigenmodes and eigenfrequenices of the system are derived and several examples 
are presented to illustrate the effect of different boundary conditions on the 
exponential decay rate of the eigenmodes. The resultant implications for the 
stability and convergence rates of finite difference computations are 
discussed. 
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while the author was in residence at the Institute for Computer Applications 
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INTRODUCTION 


In finite difference calculations of steady-state subsonic solutions of 
quasi-one-dimensional and two-dimensional Euler equations using time marching 
methods, it is often observed that when the solution has almost converged to 
steady-state the remaining residual is due to the propagation of low frequency 
waves up and down the domain. These waves are largely unaffected by numerical 
viscosity and are dissipated through the interaction with the inflow and 
outflow boundary conditions. The purpose of this paper is to examine this 
process by analyzing the unsteady linearized perturbations of a one- 
dimensional, steady, uniform, subsonic flow. For this linear problem with 
constant coefficients it is possible to derive the exact eigenmodes and 
eigenf requencies of the initial boundary value problem. This is the classical 
technique used to analyze physical and acoustical vibrations in a finite 
domain [5] and more recently used in numerical analysis to examine the 
P-stability of finite difference approximations to scalar equations [2,3]. 
The exponential decay rate of the physical eigenmodes is computed for several 
different sets of boundary conditions commonly used in finite difference 
calculations and the implications for the stability and convergence rates of 
these calculations are discussed. 

The wellposedness of both the initial boundary value problem (i.b.v.p.) 
and the steady-state boundary value problem (b.v.p.) is discussed briefly. 
The definitive analysis of the i.b.v.p. for multi-dimensional hyperbolic 
systems is given by Kreiss in [4]. Oliger and Sundstr’om [7], use an energy 
method to establish sufficient conditions for the wellposedness of the Euler 
i.b.v.p. Finally, the wellposedness of the steady-state solution to the 
nonlinear quasi-one-dimensional and two-dimensional Euler equations will be 
discussed in a forthcoming paper by Wornom and Hafez [8]. 
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2. ANALYSIS 

The equation for the unsteady linearized perturbation of a steady, uniform 
one-dimensional flow is, 



where p, u, p are the perturbation density, velocity and pressure and 
p, u, p are the steady, uniform values. 

The analysis is greatly simplified by defining the following non- 


dimensional variables 

P = P/P (2) 
u = u/c (3) 
p = p/p c 2 (4) 
x - X/L (5) 
t = Tc/L, (6) 


_ l / 0 

where c = [yp/p] z is the speed of sound. L is the physical length of the 
domain considered, so in the non-dimensional domain the subsonic inflow is 
at x = 0 and the outflow is at x = 1. 

The resultant non-dimensional equation is 
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where 

/m 1 0 

A = I 0 M 1 

\0 1 M 

and M is the Mach number of the unperturbed flow. 

Equation (7) has wave-like solutions 




= exp i(kx - u)t) U 


provided 


(kA - w!)U = 0, 


( 8 ) 


(9) 


( 10 ) 


so to/k is an eigenvalue of A and U is the corresponding eigenvector. 
The three eigenvalues of A and their corresponding eigenvectors are 



(Ha, h) 


(12a,b) 


(13a, b) 


A general eigenmode of the initial boundary value problem can be written 
as a sum of the three eigenwaves, 

-ftoi- r i(u)/A )x i(w/A„)x i(w/A„)x 

U = e [ctje U t + « 2 e U 2 + « 3 e U 3 J . (14) 
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The eigenfrequency oj and the values of the constants , a a ^ are 
determined by the three boundary conditions. 

At the inflow boundary at x = 0 there are two boundary conditions which 
when linearized and non-dimensionalized have the form. 


C in U = °> 


(15) 


where C^ n is a 2 x 3 matrix. Substitution of (14) into (15) yields the 
equation, 


where 




c i„( D i 



(16) 


(17) 


A necessary condition for the initial boundary value problem to be 
wellposed is that the 2 x 2 matrix 



is nonsingular and so can be inverted to obtain and , the values of 

the incoming characteristics, as a function of a^, the value of the outgoing 


characteristic 
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Similarly the outflow boundary condition yields one equation of the form 


C U = 0 
out 


and substitution of (14) produces 



where 


( 18 ) 


(19) 


( b 31 b 32 



i(w/X 1 ) 

e 


0 , 


i(co/X 2 ) i(u>A 3 ) 

U 2 e U 3 ). 


( 20 ) 


The second necessary condition for the wellposedness of the initial 
boundary value problem is that bgg is nonzero so that the value of the 
incoming characteristic can be determined as a function of and the 
values of the outgoing characteristics. 

Equations (16) and (19) can be written jointly as 


B(oj) | a 2 


= 0 . 


To obtain a nontrivial eigenmode B(as) must be singular 
(a^ a 2 «g) T must be a corresponding null vector, 
frequenices can be calculated from the following determinant 


( 21 ) 


and the vector 
Thus the eigen- 
equation 


det B(w) = 0. 


( 22 ) 
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The matrix B can also be used to examine whether the steady-state 
boundary value problem is wellposed. The three requirements for wellposedness 
are that a solution exists, is unique, and small perturbations in the boundary 
data produce small perturbations in the solution. 

The linearized steady-state boundary value problem has a zero solution and 
this solution is unique provided there are no nonzero solutions to 



( 23 ) 


i.e., provided that B(0) is non-singular. 

A perturbation of the boundary data leads to an equation of the form 


f; 

w v 3 

which, provided B(0) is nonsingular, can be solved to obtain 
T 

(a^ dg a g) which define the characteristic perturbations of the steady- 

state solution. 

Thus the linearized steady-state boundary value problem is wellposed if, 
and only if, det B(0) is nonzero, or alternatively the initial boundary value 
problem does not have a zero eigenf requency . 



3. EXAMPLES 

(a) Entropy, Enthalpy Specified at Inflow, Pressure at Outflow 


The physical boundary conditions are 
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X 

n 

o 

( p'/p' y - p/p 

r : 1 u' 2 + Yp ; - 

1 2 p 

Y 

Y-l -2 , 

= 2 U + 

YP 

P 

(25a) 

(25b) 


X = L 

P' = P> 



(25c) 

where p', u', p' 

are the 

unsteady physical 

variables 

which are 

a sum of the 

steady-state and 

unsteady perturbation 

variables . 

The 

corresponding 


linearized non-dimensionalized equations are 


x = 1 (0 0 1) 

At x = 0 substitution of the eigenvector definitions (lib), (12b), and (13b) 
into the eigenmode definition (14) yields 





Substitution of this equation into (26a) produces the characteristic inflow 
boundary condition 
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i\ Z 1 1 l \ /“ a 

J ° 1 ^ “2 

\° 1 l ' \ a 3! 



(y-l)(l~M )) 


0. 


Similarly at x = 1 



= e 


-iwt 


/o^ exp(ico/X 1 ) 
«2 exp(io)/X 2 ) 
exp(ia>A 3 ) 


and substitution into (26b) produces the characteristic outflow 
condition 


(0 0 1) 



ot^ exp(iw/X 3 ) 
a 2 exp(ia)/X 2 ) 
a 3 exp(iw/X 3 )/ 


- (o 


exp(iw/X 2 ) 


exp(iw/X_)) 


= 0 


Together equations (28) and (30) define the matrix B 



(y-l)(l+M) 

exp(iw/X 2 ) 


0 

(y-l)(l-M) 

exp(ico/X 3 )y 


(28) 


(29) 


boundary 


(30) 


( 31 ) 
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The eigenfrequencies are given by 


det B = (y-1) (1-M) exp(iu/X„) - (1+M) exp(io)/A-)] = 0, (32) 


==> exp ( 


==> w = 


2iw \ 1+M 


\1 - M' 


1-M 


2/ 1 - M ’ 


1 + M\ 
1-M/ 


"i logh-^7 + 


(33) 


(34) 


where n is an integer. 

Thus there is an infinite set of discrete eigenf requencies . 
to define a decay rate 

a n d = f " 

For this example 

1 - M 2 , + M-, 

\ - — — 108 • 


It is useful 


(35) 


(36) 


The amplitude of the eigenmode grows, or decays, as exp(-at), so the 
requirement for all eigenmodes to decay is > 0 for every n. In this 

example the requirement is satisfied and so any initial disturbance at 
t = 0 will decay exponentially. 


( b ) Mass Flux, Enthalpy Specified at Inflow, Pressure at Outflow 
The physical boundary conditions are 



| p'u" = p u 


(37a) 

X = 0 < 

! V u- 2 + K 

\ 2 p 

- Y “ 1 u 2 + Yp 

(37b) 



p 


X = L 

P' 5= P- 


(37c) 
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Omitting the algebraic details the resultant matrix B is 


B = 


1 + M 
(Y-l)(l+M) 
exp(ico/X 2 ) 


M - 1 
(Y-l)d-M) 
exp(iw/Ag) 


( 38 ) 


The eigenfrequencies are 


to 

n 




log 


/(I + M) [1 + M(Y - 1)]\ 
\(1 " M)[l - M(Y ~ 1) ] / 


+ 2(n + V 2 X . 


(39) 


The decay rates are 


a 

n 


1 “ M 2 ,_/(l + M) [ 1 + M(Y - 1 )]\ 
2 10g \(l - M)[l - M(y - 1)]/ 


(40) 


(c) Density, Pressure Specified at Inflow, Pressure at Outflow 


The physical boundary conditions are 



(41a) 

(41b) 


X = L p' = p . 


(41c) 


The matrix B is 


B = 


'1 

0 

.0 


(42) 


exp(iio/X2) expCiw/A^) 
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The eigenf requenices are 


w 

n 


(1 - M )mt , 


(43) 


and the decay rates are zero. 

Since one of the eigenf requenices is zero the steady-state boundary value 
problem is ill-posed, as discussed earlier. 

(d) Density, Velocity Specified at Inflow, Pressure at Outflow 
The physical boundary conditions are 



The matrix B is 


(° 

\0 exp(iw/X2) 


exp(ico/X 3 ) ; 


The eigenf requenices are 

a) = (1 - M 2 ) (n +V 2 )ir 


(44a) 

(44b) 


(44c) 


(45) 


(46) 


and the decay rates are zero. 

In this example the steady-state boundary value problem is wellposed but 
because of the zero decay rates unsteady oscillations will continue 
indefinitely without exponential growth or decay. 
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(e) Non-reflecting Boundary Conditions 

The full nonlinear non-reflecting boundary conditions specify entropy 
and the appropriate Riemann Invariant at the Inflow, and the other Riemann 
invariant at the outflow [4] 


X = 0 


p'/p' Y = p/p Y 


2 , — 2 - 
u + — =- c = u + — r c 
Y“1 Y-l 


2 _ 2 — 

X = L u' - — =- c' = u =- c 

y-l y-l 


The matrix B is 


B = 


y-l 


0 

2 


0 

0 


exp(ia)/A^) 0 -2 exp ( iw /X ^ ) , 


(47a) 

(47b) 

(47c) 


(48) 


Det B = 0 leads to a = + 00 which reflects the fact that with these 
boundary conditions the unsteady perturbations become zero after the finite 
time it takes for all three characteristic waves to cross the domain once. 


CONCLUSIONS 

The calculation of the exponential decay rates of physical eigenmodes has 
implications for the stability and convergence rates of time-marching finite 
difference computations. If the analytic problem has exponentially increasing 
eigenmodes then for sufficiently fine grid resolution a time-accurate 
numerical solution will exhibit corresponding exponentially increasing 
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eigenmodes. In a forthcoming paper, Trefethen-*- will prove that for a linear 
constant coefficient system such as this the three conditions: 

(i) Exponentially decaying physical eigenmodes, 

(ii) Dissipative interior numerical scheme, 

(iii) GKS-stable numerical boundary conditions, 
are sufficient to ensure the P-stability of a time-marching method for a 
sufficiently fine grid. P-stability was defined by Beam, Warming and Yee [2] 
and corresponds to GKS-stability with the additional requirement that none of 
the numerical eigenmodes increases exponentially. The precise definition of 
the theorem and its proof are given by Trefethen 1 , but in essence the argument 
is that condition (i) ensures that low frequency physical waves decay, while 
conditions (ii) and (iii) ensure the decay of high frequency waves, both 
physical and non-physical. 

The exponential decay rates for the physical eigenmodes also provide a 
useful lower limit on the spectral radius of the finite difference time- 
marching procedure. If a physical eigenmode decays as exp(-ot) with 
a > 0, then for a sufficiently fine grid the corresponding numerical 
eigenmode decays approximately as exp(-o n At) where n is the iteration 
number and At is the time-step. As the grid is refined with At/Ax held 

constant, At -»■ 0 and so the spectral radius is no less than 
2 

1 - oAt + 0(At ). If a=0, as in example (d), the physical eigenmodes are 
neutrally stable and so the numerical convergence rate towards steady-state is 
due solely to numerical dissipation. If this dissipation is of n th order then 
the corresponding spectral radius is 1 - 0(At n ). Non-reflecting boundary 


Trefethen, L. N. , 1983, Courant Institute, New York University, NY, personal 
communication. 
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conditions as in example (e) clearly give a much faster rate of convergence, 
but in two or three dimensions perfectly non-reflecting boundary conditions do 
not exist and in general the best that can be achieved is that there is zero 
reflection for locally plane waves propagating in a particular chosen 
direction [1]. 

It is not clear to what extent the conclusions for this model problem, 
with linearized perturbations and constant coefficients, are valid for more 
general flows such as transonic quasi-one-dimensional and two-dimensional 
flows. Nonlinear mechanisms at sonic lines and shocks are undoubtedly very 
important. However the decay to steady-state of low frequency waves will 
still depend on the physical boundary conditions and so this analysis should 
provide insight into the effect of the boundary conditions. 
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